/* Black Scholes versao sequencial 
   Trabalho 1 - Programacao Concorrente

   Henrique Teixeira Arroyo 	6419241
   Julio Polzato		7288971
   Henry D Leão			5881632
 
*/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <limits.h>

#define NSUM 25

// Function returns mean value of M numbers pointed by values
float mean(float *values, int M){
	int i;
	float sum = 0;
	
	for(i=0;i<M;i++){
		sum+=values[i];
	}
	
	return sum/(float)M;
}

// Function returns standard deviation value of M numbers pointed by values
float stddev(float *values, float mean, int M){
	int i;
	float deviation,sumsqr = 0;
	
	for(i=0;i<M;i++){
		deviation = values[i] - mean;
		sumsqr+=deviation*deviation;
	}
	
	return sqrt(sumsqr/(float)M);
}

// Myrandom implementation
const unsigned int A = 1664525, C = 1013904223;

unsigned int myrandom(unsigned int x){
        return x*A+C;
}

// Gaussian distribution random number generator
double gaussrand(unsigned int *seed){

	double x = 0;
	int i;
	
	for(i = 0; i < NSUM; i++){
		*seed = myrandom(*seed);
		x += (double)*seed/(double)UINT_MAX;
	}
	
	x -= NSUM / 2.0;
	x /= sqrt(NSUM / 12.0);

	return x;
}

// Main function of black scholes, get inputs and run the algorithm
int main(){

	float S,E,r,sigma,T,t,*trials,random,a,b,e;
	float media,desvpad,confwidth,confmin,confmax;
	int M,i;
	unsigned int seed;

	scanf("%f",&S);
	scanf("%f",&E);
	scanf("%f",&r);
	scanf("%f",&sigma);
	scanf("%f",&T);
	scanf("%d",&M);
	
	trials = (float *)malloc(sizeof(float)*M);
	
	seed = time(NULL);
	
	a = (r-(sigma*sigma/2))*T;
	b = sigma*sqrt(T);
	e = exp(-r*T);
	
	// Run trials	
	for(i=0;i<M;i++){
		
		// get trial i
		t = S*exp(a+b*gaussrand(&seed));
		trials[i] = e*fmax(t-E,0);
	}
	
	media = mean(trials,M);
	desvpad = stddev(trials,media,M);

	confwidth = 1.96*desvpad/sqrt(M);
	confmin = media - confwidth;
	confmax = media + confwidth;

	printf("S\t%f\nE\t%f\nr\t%f\nsigma\t%f\nT\t%f\nM\t%d\n",S,E,r,sigma,T,M);
	printf("Confidence interval: (%f, %f)\n",confmin,confmax);
	
	free(trials);
	
	return 0;
}
